AD-A120  096  NORTHWESTERN  UNIV  EVANSTON  IL  DEPT  OF  MECHANICAL  AND— ETC  F/G  20/9 
COMPUTER  SIMULATION  OF  KINETIC  PROPERTIES  OF  PLASMAS. (U) 

AUG  82  J  DENAVIT  .  N00014-75-C-0473 

liNfl  ASrIFIFD  9 _ _ _  NL 


AD  A120096 


SECl'R  t>  Classification  0»  TniS  pacL  f*y>en  Pete  Entered,  Unclassii ieG 


O- 

Z) 

O 

■J 

J 


!ba 


REPORT  DOCUMENTATION  PAGE 

READ  INSTRUCTIONS 

BEFORE  COMPLETING  FORM  \ 

I  1  «£  P OR T  HUMBER 

9 

2  GOV T  ACC  l  SSIDN  NO- 

3  RtTClPIENVi  CATALOG*. 

c 

*  Title  (mnd  Subtitle) 

COMPUTER  SIMULATION  OF  KINETIC  PROPERTIES 

OF  PLASMAS 

1 

S  TYPE  or  REPORT  A  PE^ 

FINAL 

sm 

6  performing  org.  report  number 

7.  AuTMORfij 

J.  Denavit 

•  CONTRACT  OR  GRANT  NUMBERS 

N00014-75-0473 

9  PERFORMING  organization  name  and  address 

Northwestern  University 

Evanston,  IL  60201 

1C  PROGRAM  ElEmEkT.  PROJECT.  TASK 

AREA  6  WORK  UN(T  NUMBERS 

11  CONTROLLING  Office  NAME  AND  AOORESS 

Office  of  Naval  Research,  Physics  Program 
Arlington,  Virginia  22217 

)2  REPORT  DATE 

August  1982 

13.  NUMBER  or  PAGES 

14  monitoring  AGENCY  NAME  A  ADDRESS^//  dllteront  from  Controlling  Ollleo) 

15.  SECURITY  CLASS,  (ot  thle  report) 

Unclassified 

Mm  DECL  ASSIFIC  ATION/ DOWN  GRADING 
SCHEDULE 

16  Distribution  statement  (oi  thi,  R&port) 

Approved  for  public  release;  distribution  unlimited 

• 

17.  DiSTPibuTion  STATEMENT  (©  t  the  mbetrmct  entered  In  Block  20,  if  different  from  Report) 

DTI  C 

^  ELECTED 

IB  SUPPLEMENTARY  NOTES  PI  ; r7!T 

'  OCT  1  2  1982 

r".  i 

A 

19  KEY  WORDS  (Continue  on  reveree  aide  If  neceeemry  end  Identify  by  block  number) 

Computer  simulation,  plasma  ion  beams,  whistler. 

/ 

-1 _ _ _ 

20..  Abstract  (Continue  on  revet ee  elde  It  neceeemry  mnd  Identity  by  block  number) 

The  research  was  directed  toward  the  development  and  testing  of  new  numerical 
methods  for  particle  and  hybrid  simulation  of  plasmas,  and  their  application 
to  physical  problems  of  current  significance  in  space  physics  and  in  energetic 
charged-particles  systems.  This  project  will  terminate  on  August  31,  1982 
and  this  Final  Report  describes  (I)  the  research  accomplished  since  the  last 
renewal  on  April  1,  1981  and  (II)  a  perspective  of  the  work  done  since  the 
beginning  of  the  project  in  February  1972. 

- - - 

S/N  0103-014-  6601 


82  10  12 


CLASSIFIED 


lCI*T§’  CLASSIFICATION  OF  THIS  P  AGE  (*»*•*  Pete  Entered) 


FINAL  REPORT 


To:  OFFICE  OF  NAVAL  RESEARCH  (Physics  Branch) 

Title :  COMPUTER  SIMULATION  OF  KINETIC  PROPERTIES  OF  PLASMAS 
Contract  No.  N00014-75-C-0473 

Principal  Investigator :  Jacques  Denavlt ,  Professor 

Department  of  Mechanical  and 
Nuclear  Engineering 
Northwestern  University 
Evanston,  Illinois  60201 


Date:  August  1982 


TABLE  OF  CONTENTS 


Page 


INTRODUCTION . 1 

RESEARCH  ACCOMPLISHED  SINCE  LAST  RENEWAL . 2 

PERSPECTIVE  OF  WORK  DONE  SINCE  BEGINNING  OF  CONTRACT . 5 

PUBLICATIONS . 12 

OTHER  RESEARCH  FUNDS . 14 

PERSONNEL  AND  TIME  ALLOTMENT . 15 

INSTITUTION  PARTICIPATION . 15 

UNEXPENDED  FUNDS  AT  PRESENT  TERMINATION  DATE . 15 


-1- 


INTRODUCTION 

The  research  was  directed  toward  the  development  and  testing  of  new 
numerical  methods  for  particle  and  hybrid  simulation  of  plasmas,  and  their 
application  to  physical  problems  of  current  significance  in  space  physics  and 
in  energetic  charged-particles  systems.  This  project  will  terminate  on  August 
31,  1982  and  this  Final  Report  describes  (I)  the  research  accomplished  since 
the  last  renewal  on  April  1,  1981  and  (II)  a  perspective  of  the  work  done 
since  the  beginning  of  the  project  in  ftebruary  1972. 
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I.  RESEARCH  ACCOMPLISHED  SINCE  LAST  RENEWAL 
During  the  present  period,  research  on  the  project  has  been  concerned 
with  (A)  particle  diffusion  due  to  stochastic  orbits,  (B)  development  of 
field-implicit  hybrid  codes;  (C)  computer  studies  of  ion  beams  and  rings;  and 
(D)  saturation  of  the  bump-on-tall  instability. 


A.  Particle  diffusion  due  to  orbit  stochastlcity. 

Particle  simulation  studies  of  diffusion  due  to  stochastic  orbits  have 
been  extended  to  general  mode  structures  with  arbitrary  electric  field 
amplitudes  and  phases.  This  general  case  requires  the  commutation  of  particle 
trajectories  by  numerical  integration  methods  (e.g.  leap-frog),  which  fail  to 
yield  convergent  values  of  the  position  and  velocity  of  a  stochastic  particle 
at  a  given  time.  This  loss  of  accuracy,  which  is  the  result  of  rapid 
amplification  of  truncation  and  round-off  errors  by  the  K-S  entropy,  cists 
some  doubts  on  particle  simulations  specifically  aimed  at  stochastic 
effects.  To  what  extent  are  stochastic  effects  observed  in  particle 
simulations  true  physics,  rather  than  consequences  of  computation  errors? 

This  difficulty  may  be  alleviated  by  using  trajectory  models,  such  as  the 
"standard  mapping",  which  do  not  require  numerical  integration.  However, 
these  models  correspond  to  particular  amplitudes  and  phases  and  have  been 
found,  in  our  earlier  computations,  to  yield  unphysical  regularities. 

A  study  of  the  convergence  of  stochastic  orbit  computations  was  therefore 
carried  out.  Since  only  time  or  aggregate  averages  are  important  in  the 
applications  considered,  the  convergence  of  numerical  integration  was  tested 
in  terms  of  these  global  quantities,  rather  than  in  terms  of  the  detailed 
motions  of  individual  particles.  Bbr  example,  a  stochastic  orbit  tends  to 
fill  an  entire  region  of  phase  space,  with  the  particle  visiting  every  part  of 
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the  region  If  sufficient  time  is  allowed.  These  regions  have  complex 
topologies  and  boundaries  which  depend  on  the  amplitudes  and  phases  of  the 
electric  field.  This  study  showed  that  sufficiently  accurate  numerical 
integrations  can  generate  such  regions  with  as  much  accuracy  as  desired, 
although  the  phase  points  of  the  particle  covering  the  region  always  have  a 
large  erratic  component.  This  research  is  presented  in  the  paper  "Convergence 
of  Stochastic  Orbit  Computations,"  by  I.  Dilber,  J.  M.  Walsh  and  J.  Denavit, 
which  has  been  submitted  for  publication  to  the  Journal  of  Computational 
Physics. 

B.  Development  of  field-implicit  hybrid  Monte-Carlo  codes. 

Recent  particle  simulation  methods  with  implicit  field  computation 
(Nason,  J.  Cbmp.  Fhys.  41,  223  (1981);  Denavit,  J.  Cbmp.  Fhys.  42,  337  (1981); 
Friedman,  Langdon  and  Cohen,  Comments  Plasma  Phys.  6_,  225  (1981))  permit  the 
use  of  large  time  steps,  At>>uip^,  which  are  not  limited  by  high-frequency 
electron  oscillations.  These  methods  allow  the  development  of  new  hybrid 
codes  in  which  relatively  dense  and  cold  plasma  components  are  treated  as 
fluids,  while  energetic  electron  or  ion  are  treated  as  particles. 

Computational  studies  based  on  hybrid  codes  of  this  type  could  include 
momentum  and  energy  exchange  mechanisms,  which  are  not  possible  in  MHD  codes 
and  impractical  in  particle  simulations  because  of  computer  limitations.  Such 
codes  are  expected  to  play  an  important  role  in  future  computational  studies 
of  beam  and  microwave  generation  devices. 

During  the  past  year,  work  has  proceeded  on  a  one-dimension*. 1  electro¬ 
static  Hybrid  Monte-Carlo  (HMC)  code  designed  to  address  computational  issues 
expected  in  the  development  of  larger  two-dimensional  electromagnetic  hybrid 
codes.  In  this  code,  ions  and  thermal  electrons  are  treated  as  separate 
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fluids  following  the  equations  of  Braginskii  and  including  electron  inertia. 
Energetic  electrons  are  treated  3*  particles  interacting  with  the  fluids 
through  the  field  and  through  collisions!  energy  exchange  and  Monte-Carlo 
scattering.  Mass  exchange  between  particle  and  fluid  electrons  depends  on  the 
local  density  and  tenperature  via  the  collision  frequency,  so  that  the 
electric  fluid  can  never  develop  supra-thermal  properties  inconsistent  with  a 
fluid  treatment.  This  experimental  code  has  allowed  the  study  of  specific 
Issues  related  to  emission  and  re-absorption  algorithms,  particle  weights, 
fluid  transport  algorithms,  zoning  and  boundary  conditions.  The  results  of 
these  studies  and  examples  of  applications  of  HMC  will  be  presented  as  part  of 
the  paper  "Time-Implicit  Hybrid  Simulations  and  Stochastic  Orbit  Computations" 
by  J.  Denavit  and  I.  Dilber,  to  be  presented  at  the  24th  tanual  Meeting  of  the 
Division  of  Plasma  Riysics  of  APS,  New  Orleans,  Nov.  1-5,  1982. 

C.  Computer  studies  of  ion  beams  and  rings. 

Research  in  association  with  the  Laboratory  of  Plasma  Studies  at  Cornell 
University  has  been  continued.  This  work  concerns  computational  studies  of 
azlmuthally-symmetrlc  ion  beams  and  field-reversed  ion  rings  with  toroidal 
magnetic  fields.  Research  in  this  area  will  be  presented  in  the  paper  "Hybrid 
Simulation  of  Ion  Beam  Transport  in  Plasma  Cannels”  by  A.  Mankofsky,  R.  N. 
Sudan  and  J.  Denavit  at  the  24th  Annual  Meeting  of  the  Division  of  Plasma 
ttiysica  of  APS,  New  Orleans,  Nov.  1-5,  1982. 

D.  Saturation  of  the  bump-on-tall  instability. 

Research  on  the  saturation  of  the  single-mode  bump-on-tail  instability 
have  been  completed.  The  collisionless  case  and  weak  colllsional  effects  have 
been  considered.  In  the  collisionless  case,  a  mode  close  to  marginal 
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2 

stability  saturates  with  a  field  amplitude  E«A  ,  where  A“2(a),-w)/w,«l  defines 

the  difference  between  the  mode  frequency  m  and  the  frequency  cu*  of  the 

-3 

marginally  stable  mode.  For  weak  collision  frequencies  v~10  y.  where  Y  is 
the  instability  growth  rate,  much  larger  saturation  amplitudes,  E«A<*,  with 
q~0.6  are  obtained.  These  results  have  been  related  to  analytical  theories  of 
the  saturation  by  O'Neil  and  by  Simon  and  Rosenbluth.  This  research  is 
presented  in  the  paper  "Simulations  of  the  Single-Mode  Bump-on-Tail 
Instability"  by  J.  Denavit .  which  has  been  submitted  for  publication  to  the 
Physics  of  Fluids. 

II.  PERSPECTIVE  OF  WORK  DONE  SINCE  BEGINNING  OF  CONTRACT 

This  section  summarizes  and  evaluates  the  activities  during  the  life  of 
the  contract,  including  specifically  (A)  overall  statistics,  (B)  a  scientific 
perspective,  and  (C)  general  observations. 

A.  Overall  statistics. 

This  contract  was  initiated  on  February  2,  1972  and  has  supported  a 
plasma  computational  effort  at  Northwestern  University  jointly  with  a  contract 
from  the  Department  of  Energy,  which  was  initiated  at  the  same  time.  Both  ONR 
and  DOE  contracts  will  terminate  on  August  31,  1982. 

During  this  period  the  ONR  contract  has  supported  the  equivalent  of 
approximately  18  months  of  time  of  the  principal  investigator,  36  months  of 
research  associate  (post-doctorate)  time,  and  48  months  of  research  assistant 
(graduate  student)  time. 

Two  Ph.D.  students  have  done  a  substantial  part  of  their  research  on  the 
contract  at  Northwestern  University.  In  addition,  the  principal  investigator 
has  been  closely  associated  with  the  research  of  three  Ph.D.  students  at 
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Cornell  University. 

A  total  of  16  research  papers  in  the  Physics  of  Fluids,  the  Journal  of 
Computational  Riyscis,  Comments  in  Plasma  Physics,  and  the  Journal  of 
Geophysical  Research  have  been  published  with  sponsorship  from  the  contract 
(see  Publication*). 

B.  Scientific  perspective. 

The  general  theme  of  the  research  under  the  contract  has  been  on  the 
development  of  new  computer  simulation  methods  and  their  application  to 
plasmas  of  interest  in  space  physics  and  in  energetic  charged  particles 
systems.  These  studies,  patlcularly  those  related  to  magnetospheric  Very-Low- 
Frequency  emissions  and  to  computer  simulations  of  ion  beams ,  were  done  in 
close  cooperation  with  the  Laboratory  of  Plasma  studies  at  Cornell  University. 

The  research  carried  out  under  the  contract  since  its  initiation  may  be 
classified  according  to  whether  the  major  contribution  was  to  (1)  computer 
simulation  methods  or  (2)  elucidation  of  new  physical  phenomena. 

1.  Computer  simulation  methods. 

a)  Noise  in  particle  simulations,  including  understanding  of  how  noise  is 
generated  (refs.  2,  15),  methods  to  reduce  noise  (and  incude 
colllsional  effects)  by  periodic  re-initialization  using  a 
reconstructed  distribution  function  (ref.  3)  and  non-random 
initializations  (ref.  16). 

b)  Stochastic  particles  in  computer  simulations,  including  the 
identification  of  stochastic  particles  in  ion  ring  particle 
simulations  and  of  stochastic  single-particle  instabilities  in 
linearized  particle  codes  (ref.  13),  and  the  study  of  convergence  of 


stochastic  orbit  computations  (ref.  18).  The  occurrence  of  stochastic 
single-particle  instabilities  is  an  unresolved  problem,  which  remains 
a  major  obstacle  to  further  development  of  linearized  particle  codes. 


c)  The  development  of  particle  simulations  with  large  time  steps,  limited 
by  the  condition  WgAt<<l,  where  u>g  is  the  bounce  (or  trapping) 


frequency,  rather  than  by  <ypAt<<l  or  u)cAt<<l,  where  u>p  is  the  plasma 


frequency,  and  <»>c  is  the  electron  gyrof requency.  Since  a»g«u)p,  ioc 
for  most  applications  these  methods  considerably  extend  the  domain  of 
application  of  particle  simulations.  Two  separate  concepts  were 
developed:  (i)  Time  steps  in  which  high-frequency  oscillations  have 
been  Integrated  and  which  advance  only  slowly  varying  amplitudes  and 
frequencies  (ref.  9).  This  approach  has  been  applied  to  studies  of 
magnetospheric  VLF  emissions  (refs.  10,  12  and  14)  to  a  study  of  the 
saturation  levels  of  the  bump-on-tail  instability  (ref.  17),  and  by  J. 
L.  Vomvoridis  and  other  authors  to  electron  cyclotron  oscillation  for 
microwave  generation,  (ii)  Time  steps  in  which  high-frequency 
oscillations  are  filtered  out  (ref.  11).  This  approach  has 
potentially  very  Important  applications.  It  was  also  developed 
independently  by  Mason  and  is  now  studied  and  extended  by  several 
other  authors. 


Elucidation  of  new  physical  phenomena. 

a)  Fundamental  studies  of  the  whistler  sideband  instability  (ref.  6)  and 
of  the  nonlinear  elongation  of  electrostatic  and  whistler  wavepackets 
due  to  trapping  of  resonant  particles  (ref.  5).  In  the  electrostatic 
case,  this  theory  of  wavepacket  propagation  has  been  the  object  of 
experimental  verification  in  laboratory  plasmas  (N.  Sato,  K.  Saeki  and 
R.  Hatakeyama,  Phys.  Rev.  Lett.  38,  1480  (1977)). 
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b)  Numerical  and  analytical  studies  of  electron  heating  due  to  parametric 
instabilities.  Including  the  effect  of  collisions  to  suppress  the 
generation  of  suprathermal  electrons  (ref.  7). 

c)  Comprehensive  theoretical  studies  of  magnetospheric  Very  low  Frequency 

(VLF)  emissions ,  including  particle  correlations  due  to  interaction 

with  a  whistler  wave  propagating  in  a  non-uniform  medium  (ref.  10), 

and  the  determination  of  the  resulting  nonlinear  growth  rates  (ref. 

12).  A  theory  of  magnetospheric  VLF  emissions  must  account  for  the 

following  features  (i)  the  triggering  of  monochromatic  emissions  by 

signals  of  sufficient  strength  and  duration,  while  the  background 

noise  and  weak  short  signals  are  not  amplified,  and  (ii)  the 

occurrence  of  frequency  changes  after  the  emissions  have  reached 

sufficiently  large  amplitude.  A  nonlinear  mechanism  exhibiting  these 

features  with  fixed  and  varying  frequencies  was  identified  and 

examined  analytically  and  by  computer  simulation  techniques.  This 

mechanism  depends  on  a  simultaneous  propagation  and  amplification  of 

wave  packets  along  geomagnetic  lines  to  maintain  the  nonuniformity 

ratio  R«VBq/Bw  in  the  regime  |r|-0.5,  corresponding  to  maximum 

amplification.  (BQ  is  the  geomagnetic  field  and  Bw  is  the  wave 

magnetic  field).  Fbr  a  constant  frequency,  this  condition  yields 

triggering  thresholds  related  to  the  properties  of  the 

magnetosphere.  Fbr  a  varying  frequency  w(t),  it  yields  the 
2 

condition  3<u/3t*u)t  for  the  large-amplitude  portion  of  the  risers, 

1/2 

where  denotes  the  trapping  frequency  of  the  wave  (ref.  14). 

C.  General  observations 


During  the  past  ten  years,  computer  simulation  of  plasmas  has  been  a 
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rapidly  growing  field  of  physics  and  has  played  an  essential  role  in  the 
development  of  plasma  theory.  In  addition  to  this  research  role,  which  is 
expected  to  continue,  computer  simulation  is  also  becoming  a  cost-effective 
design  tool  to  provide  accurate  performance  predictions  for  large  plasma 
devices,  which  are  now  entering  the  engineering  phase. 

Computer  simulation  of  plasmas  is  generally  thought  today  of  comprising 
two  distinct  areas : 

(i)  Fluid  (or  MHD)  simulation,  and 

(ii)  Kinetic  (including  particle,  Vlasov  and  Fokker-Planck) 
simulation. 

Fluid  simulation  has  generally  been  applied  to  large-scale  problems  directly 
related  to  the  behavior  of  experimental  devices.  Kinetic  simulation,  on  the 
other  hand,  has  generally  been  applied  to  basic  physical  problems  in  which  the 
particle  distributions  deviate  significantly  from  a  local  Maxwellian 
distribution,  such  as  when  wave-particle  resonances,  trapping  or  stochastic 
heating  occur.  While  kinetic  simulation  usually  cannot  be  directly  applied  to 
large-scale  problems ,  their  results  provide  phenomenological  models  which  can 
be  incorporated  to  some  extent  into  fluid  simulations. 

Unfortunately,  this  dichotomy  gives  a  rather  imperfect  representation  of 
plasmas  of  interest  in  space  applications  or  in  microwave  or  beam 
generation.  In  most  cases  these  plasmas  retain  ^trong  non-equilibrium 
properties,  which  cannot  be  incorporated  into  fluid  models,  yet  occur  on  long 
time  and  space  scales,  beyond  the  reach  of  conventional  kinetic  simulation. 
Fhture  progress  in  computer  simulation  of  plasmas  ^11  therefore  depend  to  a 
significant  extent  on  the  development  of  new  hybrid  codes,  in  which  dense  and 
relatively  cold  plasma  components  are  treated  as  fluids,  while  energetic 
components  are  treated  kinetlcally  (e.g.  as  particles).  Such  codes  are 
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necessary  to  treat,  for  example,  Interactions  between  electron  or  ion  beams 
(or  rings)  and  plasmas,  suprathermal  components  or  wall-plasma  interaction 
phenomena.  Since  a  computation  can  operate  only  at  a  single  time  scale,  high- 
frequency  phenomena  (such  as  plasma  oscillations)  must  be  eliminated  from 
these  hybrid  simulations.  The  recent  time-filtering,  field-implicit  particle 
simulation  methods  are  therefore  an  important  step  in  the  development  of  such 
codes. 

Another  area  of  future  development  in  particle  simulations  concerns 
linearised  particle  codes.  Most  plasma  devices  (or  experiments)  are  based  on 
two-dimensional  equilibria,  which  are  potentially  subject  to  three-dimensional 
instabilities.  It  appears  difficult  in  the  forseeable  future  to  perform  fully 
three-dimensional  particle  simulations  on  large  space-time  scales,  but 
linearized  particle  codes,  which  compute  three-dimensional  perturbations  of  a 
two-dimensional  configuration  are  practical.  Such  codes  yield  at  least  the 
initial  Instability,  and  in  some  cases  can  model  saturation  resulting  from 
phase  mixing  of  betatron  oscillations  for  example.  However,  stochastic 
instability  of  single  particle  orbits  has  been  an  obstacle  to  the  development 
of  such  codes.  In  essence  linearized  particle  codes  follow  the  separation 
between  two  states  of  a  particle  (perturbed  -  unperturbed)  and  this  separation 
grows  exponentially  for  stochastic  orbits,  rapidly  masking  collective 
phenomena  of  interest.  A  method  to  stabilize  such  single-particle  modes  would 
therefore  be  an  important  advance. 

In  the  area  of  magnetospheric  VLF  emissions  the  greatest  need  appears  to 
be  the  acquisition  of  sufficiently  reliable  data  against  which  the  relevance 
of  theoretical  models  could  be  tested.  Of  particular  importance  would  be  the 
simultaneous  in  situ  measurement  of  the  cold  electron  density,  the  resonant 
electron  distribution  (in  Che  10-30  keV  range),  and  the  pulse  amplitude  just 
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before  the  growth  In  time  occurs.  Aaong  theoretical  issues  of  interest  in 
magnetospherlc  VLF  emissions  which  remain  unresolved,  is  an  explanation  for 
the  initial  generation  of  variable-frequency  waves.  Our  theory  provides  a 
model  for  amplification  of  such  waves,  whicu  allows  predictions  of  their 
spectral  shape  in  terms  of  magnetos pheric  properties.  However,  the  initial 
generation  of  such  waves  appears  to  require  a  sideband  instability  process 
taking  place  when  an  on-frequency  emission  has  reached  a  sufficiently  large 
amplitude. 
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Research  carried  out  since  the  present  Navy  contract  started  has  resulted 
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